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Abstract 

We present a method to calculate to very high precision the coefficients of the diver- 
gences occuring in two-loop diagrams for a massive scalar field on the lattice. The ap- 
proach is based on coordinate space techniques and extensive use of the precisely known 
Green's function. 
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1 Introduction 



Lattice regularization is a convenient tool to study both the perturbative and non-perturbative 
aspects of the underlying theory. Within this framework the action is formulated on a dis- 
cretized space-time lattice and the partition function can be evaluated numerically, e.g., by 
employing Monte Carlo simulations. This feature of lattice regularization allows one to study 
the non-perturbative regime of a theory and is thus of particular interest for QCD at low 
energies. 

Eventually, however, one is interested in performing the continuum limit by letting the 
lattice spacing approach zero. In order to establish the correct link between numerical lattice 
simulations and the physical continuum, perturbative lattice calculations prove useful. They 
are, e.g., necessary to calculate the renormalization of the couplings in the Lagrangian or 
the renormalization factors of operator matrix elements. In general, it allows one to study 
perturbatively how simulation results approach continuum results. The discretization effects, 
i.e., the corrections due to the finite lattice spacing, can also be investigated and this knowledge 
may help to reduce the pertinent systematic error in lattice results. 

The explicit calculation of loop diagrams in lattice regularization turns out to be more 
tedious than in conventional continuum regularization schemes such as dimensional regulariza- 
tion. In this respect, coordinate space methods have proven useful in the evaluation of Feynman 
diagrams and allow a very precise determination of the continuum limit of one- and two-loop 
integrals [1]. In [1] this technique was applied to massless propagators. In [2] the formalism 
of [1] was extended to massive scalar fields by presenting an efficient method to calculate the 
associated Green's function to very high precision. 

In the present investigation, we will show that the precise knowledge of the Green's function 
is helpful in extracting the divergent pieces of two-loop diagrams for massive scalar fields to 
very high accuracy. The method developed is therefore suited to calculate the renormalization 
of the bare Lagrangian parameters at the two-loop level. 

This work is organized as follows. First, external momenta in the loop propagators are 
eliminated by applying the BPHZ scheme on the lattice so that the discussion can be reduced 
to sunset-type loops without external momenta in the propagators. We illustrate the BPHZ 
method in the next section for the basic sunset diagram. Then, the divergences of the reduced 
integrals are evaluated in Sec. 3 by applying coordinate space methods. In Sec. 4 we present 
the numerical results for the quadratic and logarithmic would-be divergences of the sunset 
diagram and in Sec. 5 we provide additional relations between the logarithmic pieces and one- 
loop tadpole diagrams. We summarize our findings in Sec. 6. Technicalities are deferred to the 
appendices. 

2 BPHZ scheme on the lattice 

Within the BPHZ scheme the ultraviolet divergences of integrals are removed by performing 
subtractions directly in the integrand of a Feynman integral [3]. This is achieved by subtracting 
the first few terms of the Taylor expansion in the external momenta. Originally this method was 
developed for the continuum formulation of field theories, but afterwards extended to lattice 
field theories as well [4]. 

In the present work we are particularly interested in isolating the divergent pieces of a 
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Figure 1: Shown is the sunset diagram with three massive propagators. 

given lattice loop integral, while the remaining finite pieces can be calculated numerically for 
a given lattice spacing. We restrict ourselves to massive scalar field theories on the lattice. 
For illustration the method will be worked out explicitly for the basic sunset diagram which 
arises in 4 theory at the two- loop level, see Fig. 1, but a general two- loop diagram can be 
evaluated in the same manner. The first step in the approach is the application of the BPHZ 
scheme so that the would-be divergences which occur as the lattice spacing goes to zero are 
contained in lattice integrals without external momenta, while the momentum-dependent but 
finite remainder can be neglected for renormalization purposes. Provided the lattice theory 
satisfies the power counting conditions given in [5], which is the case here, the continuum limit 
of these momentum-dependent finite pieces exists, is independent of the details of the lattice 
action and is identical with the BPHZ finite parts obtained from the corresponding continuum 
action [4]. In order to illustrate the method, we will confirm explicitly for the case of the sunset 
diagram that this remainder is indeed finite in the continuum limit. For definiteness we work 
with a hypercubical lattice A of infinite volume and the standard free scalar propagator, while 
the results are presented in lattice units if not otherwise stated. 

The basic sunset diagram for a massive scalar field with external momentum p is given by 
Fig. 1, 

r d 4 k d 4 q r d 4 k d 4 q 1 

J_„ (2tt) 4 (2tt) 4 Si > q]P) ~ J w (2vr) 4 (2tt) 4 A(k) A(q) A(p + k + q) ' U 

where the subscript S of I denotes the full sunset diagram and A is the inverse lattice propagator 

4 4 

A(p) = m 2 + p 2 = m 2 + = m 2 + 4^sin 2 (^j . (2) 

/i=i /i=i 

Let 71,72,73 be the three subdiagrams of S, Fig. 2, then the corresponding set of forests 
is given by J-'(S) = {0, S, 71, 72, 73, Sji, Sj2, £73} with the superficial degrees of divergence 
d(S) = 2 and ^(7) = otherwise. The integrand Is of the sunset diagram leads to would-be 
divergences in the continuum limit as can be seen by explicitly inserting the lattice spacing in 
Eq. (1). A convenient renormalization procedure to handle these divergences is provided by 
the forest formula of Zimmermann [3] wherein a subtraction operation is directly applied to 

Figure 2: Subdiagrams 71,2,3 of the sunset diagram. 
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the integrand. Within this approach the renormalized Feynman integrand Rs is given by 

R s = H IlHfVs, ( 3 ) 

where are lattice subtraction operators of order £#(7) as introduced in [4]. In the continuum 
limit the operators coincide with the first ^(7) terms of the Taylor expansion of the 
integrand J 7 in the external momentum of the subgraph. For the sunset diagram one obtains 

Rs = (i-4-ti- 1% 2 - ti + e s + 4 + 4] ) i s 
= (i-4-[i-4]K + 4+<J)^- ( 4 ) 

Since [1 - t|] (t^ + 1° 2 + t° l3 )I s = 0, this leads to 

R S =[1- 4] / s . (5) 
A convenient choice for the subtraction operator t s is given by [4] 

d d . 



p=0 



with 

4 4 



(6) 



t 2 s Is(k,q;p) = l s (k,q;0) + ^sinp M (JLlg\ + I sin sin p^ ( 
where the partial derivatives are to be taken at p — 0. This yields 

+4^sin 2 p^ sin 2 (/c + g)^A- 2 (A; + g)) . (7) 
Utilizing hypercubical symmetry of the lattice it can be rewritten as 

/ 1 77?^ 

4 /*(*, g; p) = A-\k)A-\q)A-\k + q)[l + -p 2 - (1 + — )fA~ 1 (k + q) 

+4E?ira; A_2 ( fc +?)) ( 8 ) 



p 2 = Y,pI = E sin2 ^- (9) 

/1=1 /1=1 

Neglecting pieces that do not lead to divergencies in the continuum limit this simplifies to 

t| I s (k, q- p) = A-\k)A-\q)A-\k + q)(l- p 2 A~\k + q) 

+AY,Pl{kTq)l^\k + q)) + ... . (10) 
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According to [4] integration of Rs in Eq. (5) yields a finite continuum limit which is independent 
of the details of the lattice action. In order to demonstrate the finiteness of the renormalized 
integral, it is thus sufficient to consider the corresponding integral in the continuum formulation 



R c s(k,q; P ) 



(2tt)* (2tt)* 
r°° d 4 k d 4 q 



D- 1 (k)D- 1 (q)- 1 D- 1 (p + k + q) 



-D-\k)D-\q)D-\k + q) (l - p 2 D" 1 {k + q) + 4 + q) 2 l D~ 2 (k + q)) j (11) 

with the inverse continuum propagator 

4 

D(p) = m 2 +p 2 = m 2 + ^2pI • (12) 

We employ dimensional regularization for the evaluation of the Feynman integral (11) which 
can then be simplified according to 



' d d k d d q „ ct 



(2ir) d (27r) d 

= £(S(S{ D " 1WB " ,(9) " lD " 1(p+ ' ;+<!) 

-D-\k)D-\q)D~ l {k + q)(l- p 2 m 2 D- 2 {k + q) + ^p 2 D^ 1 (k + ?))! + 



(13) 



where e = 4 — 0? and we have neglected finite terms. One verifies in a straightforward manner 
that the divergent pieces of the single terms in Eq. (13) indeed cancel. This is most conveniently 
done by noting that in dimensional regularization the divergent pieces of the sunset integral 
with two different masses, 

t ( 2 2 71 /r2 f d k d 1 , . 

l{m ,m ,M ;p ) = j ^ ( 27r )d( 27r )<i[ m 2 + jfc 2][ m 2 + 9 2][ M 2 + ( jfc + ? + p )2] ' 
are given by 

I dlv (rn 2 ,m 2 ,M 2 ;p 2 ) = -^- e [2m 2 (- + 1 -2 lnm 2 ) +M 2 (- + 1-2 InM 2 ) + ^1 . (15) 

The divergences of the scalar integrals in Eq. (13) can now be deduced directly from the 
expression in Eq. (15) and its differentiation with respect to M 2 . One obtains 

00 d d k d d q r . /7 . , x 

(0 , d(0 l d Rs^q-p) =0. (16) 

We have thus confirmed explicitly that the would-be divergences of the original lattice integral 
(1) are contained in the expression t 2 s Is of Eq. (10). As there is no external momentum left 
in the propagators of t$Is, one can restrict oneself in the following to loop integrals without 
external momenta in the propagators. 
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3 Coordinate space method 



After having applied the BPHZ formalism to the original sunset integral we now continue with 
the evaluation of the three lattice integrals appearing in J t 2 s Is, Eq. (10). Instead of calculating 
the divergent pieces in momentum space, it is more convenient to work in coordinate space. 
To this end, one can employ the results for the massive Green's function presented in [2]. We 
start with the integral 



d 4 k d 4 q 



w (2ir) 4 (2n) 4 A(k)A(q)A(k + q) 
where G is the lattice Green's function 



E G ^) 3 ' ( 17 ) 



I 7 (18) 



i k-x 

In lattice units the sum over G 3 is finite and can be carried out numerically (with increasing 
precision for larger values of m) utilizing the precise numerical knowledge of the Green's function 
[2]. In order to isolate the would-be divergences for a — > 0, however, one recalls that in physical 
units the sunset diagram in Eq. (1) enters with a prefactor 1/a 2 , where a is the lattice spacing 
and the physical mass \i is given by m = a/i. For large enough x the first term in the small 
mass expansion of G, i.e. the massless Green's function G m=0 , behaves asymptotically as 1/x 2 
and, hence, G z m=Q as 1/x 6 [1]. One can thus perform the sum explicitly. But for the logarithmic 
divergences which arise in the subleading mass term one observes an asymptotic 1/x 4 behavior 
of the expansion coefficient of G 3 and the summation over all lattice sites cannot be performed 
numerically. We therefore proceed by splitting the sum of the Green's functions as follows: 

J2G(xf = G(0f + £(G(s) 3 - G^(xf) + £Gi s )(*) 3 . (19) 

x<=A x+Q x+Q 

In this expression, we have introduced the continuum Green's function 

G^{x) = ^ I - 2 m\x\K 1 {m\A) , (20) 
where \x\ — Vx 2 and K\ represents the modified Bessel function of the second kind 



Kv{2) = 2 sin(un) ' = ^ i\T(y + i + Tj ' (21) 

The Green's function at the origin, G(0), has been separated from the sum, as G^(x) diverges 
at the origin. Note that we have merely subtracted (and added) the continuum Green's function 
which is sufficient to expand the difference in the first sum on the right side of Eq. (19) up 
to order m 2 . In order to be able to expand to higher orders in m 2 one must subtract further 
asymptotic pieces of the lattice Green's function which are suppressed by increasing powers in 
the lattice spacing and represent corrections to the continuum expression. They have already 
been presented in [6] , where the asymptotic behavior of the lattice Green's function for large x 
has been studied. For completeness, we present in Appendix A an alternative and much simpler 
derivation of the asymptotic form of the lattice Green's function. 
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For large x and m = the difference G(x) 3 — G^aJ(x) 3 vanishes as 1/x 8 so that the small 
mass expansion can be performed up to m 2 . To this end, we note that the small mass expansion 
of the Green's function is given by [2] 

oo oo 

G(x) = <x)m 2i + m 2 ln(m 2 ) b t {x)m 2 \ (22) 

i=0 i=0 

which converges absolutely for m < 2. The functions ai(x) can be calculated to high accuracy 
using a set of recursion relations, while the bi(x) can be given in closed form 

^ = m^^-^ b ^ x) = 3ok-^ x2)2 - 4x2+ ^ etc - (23) 

Moreover, the asymptotic large-x behavior of the coefficients a iy a® s , can be derived by expand- 
ing GaJ in m 



a?(x) 



16tt 2 



l + 2( 7s + ln(M))] + 0(l/x 2 ) (24) 



with 7^ = 0.5772 . . . being the Euler-Mascheroni constant. (In Eq. (24), subleading orders in 
1/x 2 can be deduced from the higher order corrections of the asymptotic Green's function.) 
We therefore obtain up-to-and-including order m 2 the expressions 

G(0) 3 = a (0) 3 + 3m 2 a (0) 2 ai(0) + 3m 2 ln(m 2 )a (0) 2 6 (0) + C(m 4 ) (25) 

and 

J2(G(xf - G^) 3 ) = $>o(:r) 3 - a a s (xf) + 3m 2 5> (x)V(*) - a a s {xf a? (x)) 

x^O x^O x^O 

+3m 2 ln(m 2 ) J2( a o(x)X(x) - a a s (x) 2 b (x)) + 0(m 4 ) 



= E ( a °(*) 3 - 64^) 

x^O 

+3m2 E («o(*)V*) -aSv[»» + 2 In (M) - l]) 

x^O v 7 

+3m 2 ln(m 2 ) £ (a (x) 2 6 (x) - ^JL-j) + 0(m 4 ) . (26) 

x^O 

In the continuum limit, the first and third sum in Eq. (26) lead to quadratic and logarithmic 
divergences, respectively, whereas the second sum remains finite. 

We can now turn our attention to the evaluation of the last term in Eq. (19) 

E^(-) 3 = E (ii^w^(-M)) 3 (27) 

x^O x^O 



7 



which we rewrite as 



x^O x^O 

Here, A(m, x,d) is given by 



(^-^^m^K^m^y - A(m,x,A) +J^A(m,x,4). (28) 



A(m, d) = E] C I 3 



j=0 







(• 


[3-'] 



7=0 ^ ■ 



— j,mx 



— 1 — j, m x 



3=0 











(■ 


[H 







— j,mx 



with 



C(i/,x) 



.7; 



(29) 



(30) 



and the coefficients Kq, . . . , Kg, which also depend on m and d, are given in Appendix B. The 
function A has been chosen in such a way that the difference 



b = lim 

x^O 



1 



Aix 2 x- 



-m\x 



\ 3 

Ki(m\x\)\ -A(m,x,A) 



= 0{m A ) 



(31) 



remains finite at the origin. We can thus replace the first sum in Eq. (28) by the sum over all 
lattice sites including x — 0. We are left with the evaluation of the difference on the right-hand 
side of Eq. (28) which is rewritten as 



E (iJ^ m W^i(™M)) -A(m,x,A) 
xeA L 



- b . 



(32) 



Taking in Eq. (32) the Poisson resummation 

ro ° d A p 

xeA 

yields 



/OG j4 



(33) 



E G ^) 3 = Efc? 4 

x^O 



I{m\ m% m 2 ; (2tt£;) 2 ) - A(m, 2nk, d) 



-b 



x^O 



m 



(34) 



with I the continuum sunset integral in <i-dimensional space defined in Eq. (14) and with 
external momentum 2nk, while A is given by 



A(m,p, d) = 



(47T) 



d/2 



3 m 6+2j-3d 2 ( d 

E^ r ( 3 + J - ^7^X^5+5=5 + E^T 2 + j - - 



m 



A+2j-2d 



(p2 + m 2)3+j 

m 



2 ^2 +7 •^ + E /^ i+ 7^ ( 1+ ^ 

(p 2 + m 2 ) 2+J 2 ^ 



j=0 
2+2jr'-(i 



( p 2 + m 2)l+j 



(35) 
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We have analytically continued the expressions in the first sum of Eq. (34) to rf-dimensional 
space. For d — 4, the single terms in the difference of Eq. (28) diverge at the origin and 
cannot be calculated separately, whereas for d small enough both terms remain finite and 
can be evaluated. Taking afterwards the d — > 4 limit one confirms that the divergences in d 
cancel out in the difference. By continuing to d dimensions we have automatically introduced 
dimensionally regularized integrals. In the course of evaluating sunset diagrams with momenta 
2nk in dimensional regularization it is more convenient to discuss the case with momentum 
k = separately. The results for external momentum k — are presented, e.g., in [7,8] 



I(m 2 , m 2 , m 2 ; 0) = 



3 rri 2 



1 



+ 



3 m 2 



128 7r 4 (d- 4) 2 256 7T 



(3 - 2 7E - 2 ln(m 2 ) + 2 1h(4tt)) 



1 



d-A 



+ 



m 



256tt 4 



- y - 3t| - \ - 9 1h(4tt) - 3 ln 2 (47r) + 2V3C1 2 (|) 



+ 7^(9 + 6 1h(4tt)) + ln(m 2 )(9 -67^ + 6 1h(4tt)) - 3 ln 2 (m 2 ) 



(36) 



with Cl 2 being the Clausen function. For k 7^ 0, on the other hand, one can expand in powers 
of m 2 /(2vrA;) 2 [7,8] 



Til 1 11 

l[m ,m ,m ;p 



3 m 2 



+ 



128 7T 4 (d-A) 2 



V 



+ 



12m 2 ln(m 2 ) + 6m 2 (3-27 S + 21n(47r)) 1 



5127T 4 



2048tt 4 



m 



1024vr 4 



30 



d-A 
12 lE {-3 + lE ) 



13 + A^ E - A ln(47r) + 4 ln(p 2 
+ 7T 2 + 36 1h(4tt) + 6 [ ln 2 (m 2 ) - \n 2 (p 2 ) + 2 ln(m 2 ) (-3 + 2 7e + ln(p 2 ) 
- 21n(47r)) + 21n(47r)(-27 S + ln(47r))]j + m 2 0(m 2 /p 2 ) . (37) 

For our purposes we only need to expand up to m 2 . Note that in the difference I — A the 
l/(d — A) divergences cancel and the difference I — A scales like 1/ k 6 for large k so that the 
sum over k converges. 

We are left with the evaluation of the sum A(m,x,A) in Eq. (34). To this aim, it is 

convenient to work with the integral representation of the Bessel functions K leading to, cf. 
Eq. (30), 



C(v,mx) = 



2 m 



2;/ 



d\\ 



-u-1 



exp ( - m 2 X - . 



(38) 



Since the function A is a linear combination of the C functions, it is sufficient to demonstrate 
the calculation for the sum 

Y,C{v,mx), (39) 

or in more general — if we are interested in the calculation of m 4 and higher order corrections — 
terms of the form 

^2h D (x)C(v,mx) . (40) 

Here, hn(x) is a harmonic homogeneous polynomial of degree D with the leading orders given 
by 



h (x) 



hAx] = 2x — (x 



2\2 



h e (x) = 16x b - 20xV + 5(x 2 Y, 



(41) 
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Following [1] we introduce the heat kernel 

k(X, h D ) = ^ h D (x)e" K 

X 

and rewrite the original sum in the form 



\x 2 



(42) 



1 / Air \ V f°° 2 

y2h D (x)C(u,mx) = - — / d\\-^ l e-^ x [k{l/\,h D )-h D {Q)}. 



(43) 



Poisson resummation and the fact that h D (x) is a harmonic homogeneous polynomial of degree 
D lead to the relation 

k(X, h D ) = {-l) D ' 2 \- D - 2 k{l/\ h D ). (44) 
Employing this relation we can separate the regular and singular parts of the integral 



hD{x)C{v,mx) 

x^O 



1 /4ttV f°° 

- 2 {- 2 ) I d\[k(\,h D )-h D 



v-l 



/ 2 \ 

/ m \ 



\47T/ 



2 \ 

m \ 



An J 



with the exponential integral functions E n , 

/oo 
d\e~ zX \- n . 



(45) 



(46) 



The first term in Eq. (45) is finite for all values of v and can be Taylor expanded in m 2 , while the 
exponential integral functions entail the logarithmic terms. The sum over the A functions can 
now be represented as a linear combination of regularized zeta functions and their derivatives 
with respect to s 



Us, h D ) = I d\ [y- 1 + {-i) D ' 2 \ D - s+1 ] [k(X, h) - mo)] , &\s, h D ) = ^-Us, h D ). 



The regularized zeta functions are related to the generalized zeta functions via 



x^O 



X 



2\-s 



2M0) 
r(s) Ls(s-2) 



+ £r{s, h D ) 



(47) 



(48) 



Altogether, we obtain 



x^O 



1 Wo) 
3 2 



+ m 



11 3 7£ 3 7 £ . n l 3 



2 t! 

4 ' 2 



^ + - + ^ln(47r) 



+ 3 7s 1h(4tt) - ^ ln^(47r) + 3 7e ^(2, M - 3 Ih(4tt)^(2, fc ) - 3^(2, fc ) 



+ m 2 ln(m 2 ) 
+ C(m 4 ). 



(4tt) z 



--- 3 lE + 3 ln(4vr) + 3^(2, /i ) 



— m In (m ) 



3 1 



2 (4vr) 4 
(49) 
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We can now summarize the quadratic, logarithmic and bilogarithmic divergences as well as 
the finite part of the integral in Eq. (17). From Eqs. (25), (26),(36),(37)and (49) we obtain the 
quadratic would-be divergences 



d 4 k d 4 q 



1 



quadr.div. 



a » ( °) 3+ 5:( a ° w3 -64^) + (SF 



x^O 



1 Wo) 
3 2 



(50) 



the logarithmic divergences 
n d 4 k d 4 q 



w (2tt) 4 (2tt) 4 A(k)A(q)A(k + q) liog.div 
= m 2 ln(m 2 ) 



3a (0) 2 6 (0) + ^ (3a (x) 2 b (x) 

x^O 



256vr 6 a; 4 



+ 



(47T) 



37. 



+31n(47r) + 3^(2,/i ) 

the bilogarithmic divergences 
d 4 A; d 4 q 



bilog.div. 



_ w (27r) 4 (27r) 4 A(A;)A(g)A(A; + g) 

as well as the finite remainder in the continuum limit 
d 4 k d 4 q 1 



2(4tt)' 



-m 2 ln 2 (m 2 ) , 



(51) 



(52) 



= m 



n {^Y (27r)*A(k)A(q)A(k + q) 

2 



finite 



-(2 7 £ + 21n(a;/2) - 1)) + ' 



3a (()) 2 a,(())-^(3« () (,-) 2 « l (.r) v - . , v 

37s 3 7 2 7T 2 . „ /vr\ 3, s . „ , , 3 



15 
T 



+ ^ + 2>/3Cl 2 (|) + - 1h(4tt) + 3 7s 1h(4tt) - - ln 2 (47r) + 3^(2, h ) 



-31n(47r)^(2,/i )-34 1) (2,V 



(53) 



So far, we have evaluated the divergent components of the sunset integral given in Eq. (17). 
The calculation of the divergences in the second integrand in t 2 s Is, Eq. (10), is immediately- 
obtained via 



* d 4 k d 4 q 



1 



1 d r d 4 k d 4 q 



n (2tt) 4 (2tt) 4 A(k)A(q)A(k + q) 2 3 dm 2 J_ w (2tt) 4 (2tt) 4 A(k)A(q)A(k + q) 



1 



(54) 



and by noting that this time we only need the small mass expansion at leading order due to 
the presence of the prefactor p 2 in Eq. (10). Thus, the second term in t 2 s Is contributes only to 
the logarithmic divergences. 

The third and final contribution in t 2 s Is, Eq. (10), is given by 

^^ A^tf = -iES.WEP, + VJ G , W ■ (55) 
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Here, we have introduced the notation 

-,wrwm (56) 

and the backward and forward lattice derivatives are defined by 

V*/(x) = f(x) - f(x - /i) , V„/(x) = /(x + /i) - /(x) . (57) 

One easily verifies that the most general lattice integral with one mass arising in the BPHZ 
subtraction procedure has the form 

J_ n (2tt) 4 (2tt) 4 A(k)«A(q)PA(k + q)i ' 1 j 

From Eq. (55) it is obvious how a generic integral of the form (58) can be treated with the 
proposed coordinate space method. The presented approach is therefore suited to evaluate 
any two-loop integral without external momentum which arises in the BPHZ procedure. The 
extension of this method to scalar field theories with different masses is also straightforward. 

We proceed as above by separating in Eq. (55) the contribution at the origin from the 
remaining sum and isolating the corresponding continuum result 

-\ E G ^ 2 E^ + K] 2 °^) = -\ G ^f + v ;] 2G 3(o) - [\ G ^ 2 

x ii n x^O 

x E( y M + v;] 2 g 3 (x) - gi(x) 2 (^) 2 g§(*)] - E G ^) 2 E (^r) 2 ( 59 ) 



with 

/oo t4l p ik-x 

(Note that in our notation G\ = G^J ■) The first term on the right side of Eq. (59) can be 
directly derived from the Green's function at the origin [2], and the subsequent term can also 
be calculated numerically, since the difference of functions in the sum behaves like 1/x 6 for large 
x in the m — > limit. In contrast, the rightmost term in Eq. (59) must be further decomposed 

2 



x^O /i ^ / 

= E \ G ^ x ) 2 E 2 G 3» " S(m, x, 4)] 



- lim (Gl(x) 2 <^(x) - S(m, x, 4)) + £ 5(m, x, 4) 

= J2 J rf 4 xe^[^(x) 2 ^(Ay G e (;r) _ 5(m)a;?4) j 



2 



Bm(G5(x) 2 X;(^-) (61) 



3^0 
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where the function B(m, x, d) is defined as 



B(m,x,d) =J2 X i C ( 3 



3=0 



d 



- 1 



- ./■ in .c ] + A i+4 C ( 2 
i=o 



d 



- 1 



— j,mx 



+ Y^\+iC - l-j,mx^J 



(62) 



with coefficients A.,-, which also depend on m and <i, given in Appendix B. The function 
B(m, x, d) has been introduced in order to apply the Poisson summation formula. The difference 
at the origin x = in Eq. (61) is given by 



lim(^(:r) 2 W^U G\{x) - B(m, x, 4)) = 0(m 4 



For the evaluation of the first term in Eq. (61), on the other hand, we utilize 



J d d xe ip x Gl(x) 2 J2 (J^r) ° C ^ X ) = J ^xe^GKxf 



d d 



d 2 



d d 



+ 



X 



A dm 2 8 



Gt(x) 



+ 



d-1 d 



12 dm 2 8 V<9|p| 2 \p\ d\p 



(63) 



where we introduced analytical continuation to d dimensions in order to evaluate both terms in 
the first sum of Eq. (61) separately. For k ^ we can directly apply the differential operators 
of the last equation to the result of Eq. (37). For k = 0, on the other hand, we apply the 
differential operators to the small p expansion of the continuum sunset diagram. The leading 
order of this expansion has already been presented in Eq. (36), while the p 2 terms are given by 



p 



512(-4 + rf)7r z 



+ 



-27 + 367^ + 36 ln(m 2 ) - 36 1h(4tt) + 32^012 (f ) 
18432vr 4 



(64) 



Note that for our considerations we can neglect p 4 and higher contributions. The Poisson 
resummation for B follows immediately by noting that B has the same structure as A in 
Eq. (29) and replacing re, by A, in Eq. (35). As final result we obtain 



J d 4 x [Gt(x) 2 ^(^yG C 3^)-B(m,xA)] 



9 

4 3^ 



2 Cl >(i 



(65) 



Summarizing the divergent and finite portions of the integral in Eq. (55) which remain in 
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the continuum limit one obtains 



d k d q sin (k^ + q^) 



log.div. 



(27r) 4 (27r) 4 A(k)A(q)A(k + q) 3 



2567r 6 a; 4 



+ 



^-L i (2 + 7B -ln(47r)-^(2,/i ))}, 
" d A k d 4 q sin 2 (^ + g M ) 



, (27r) 4 (27r) 4 A(fc)A(g)A(A; + ^ 



bilog.div. 



2(4tt) z 



m 2 (m ) , 



(66) 
(67) 



* d 4 k d A q sm 2 (k^ + q^) 



, (27r) 4 (27r) 4 A(A;)A(g)A(A; + g)3 



= — ao 



(0) 2 [^ + 2(a 2 (2,0,0,0)-a 2 (0)) 



^ [§ + 2 7 , + | - £ - ^ - 2 ln(4.) - 1E ln(4.) + I ln 2 (4.) - 



3£fl(2, M 



(4vr) 4 L2 ' ' 2 12 3^ 
-7i5&(2, /io) + ln(47r)^(2, /i ) + ^ 1} (2, M 

4 Numerical results 



(68) 



The Equations (50-53) and (66-68) can in principle be evaluated numerically. However, the 
respective sums converge slowly. In order to increase the precision of the numerical results one 
must therefore accelerate the convergence of the sums. This can be achieved, e.g., by using the 
zeta function technique introduced in [1]. 

The small mass expansions of the integrals in Eqs. (1) and (55) are given by 



/ 7r d k d q 1 . . 

j2^j2^ A(k)A(q)A(k + q) = g + m2 ln(m2) g + m2 ln2(m2) g 6 

oo 



m 



2i 



IT J4 



drk d 4 q sm 2 (k^ + g M ) 
n (2n) 4 (2n) 4 A(k)A(q)A(k + qy 



E h™ 21 + ln(^ 2 ) E W + ln 2 (m 2 ) E K 



in 



2i 



i=Q 



i=0 



i=0 



Recall that the expansion for the integral in Eq. (54) follows directly from the first equation 
above by differentiation with respect to m 2 . Utilizing the zeta function method we have calcu- 
lated the first few coefficients of these series. We obtain 

c = 0.0040430548122 . . . , Cl = -0.0024114634124 . . . , 

d = 0.0006968046967 . . . , e = -■ f 



and 



2 (4vr) 4 ' 



fo = 0.0000731523655 ... , # = -0.0001721159925 ... , h 



1 1 

2(4^f 



(69) 
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In principle, we can calculate higher coefficients with this method, but for renormalization 
purposes this is not necessary. Also, the precision of the coefficients can easily be improved by 
including higher orders in the zeta function technique and increasing computing time a bit. 



5 Additional relations for logarithmic terms 

The logarithmic terms can, in fact, be related to one-loop tadpoles. Starting from the identity 

r * d 4 k 1 



one has for p = 



E G (-) 2 



) 4 A(k)A(p + k) 
d 



dm 2 



G(0). 



On the other hand, decomposition of the sum leads to 

E G{xf = a (0) 2 + («o(*) 2 - + E G «(*) 2 + (™ 2 )- 



(70) 



(71) 



(72) 



x^O 

Employing the small mass expansion of the tadpole 

d 



dm 2 



G(0) = ai(0) + 6o(0) + 6o(0) ln(m 2 ) + 0(m 2 ) 



(73) 



one gets 



E ( a °(*) 2 - 16^4) = - fl o(°) 2 " E G «(x) 2 - ai(0) - 6o(0) - 6 (0) ln(m 2 ) + 0(m 2 ). (74) 

The sum over the asymptotic Green's function on the right-hand side can be calculated with 
the method presented in Sec. 3: 

E G «(x) 2 = j^~ 2 [-\-1e- ln(m 2 ) + ln(4vr) + fr(2, ho)] + 0(m 2 ). 

x^O \ ' 

Comparison with Eq. (51) leads immediately to the logarithmic divergence 



(75) 



7T JA 



d A k d A q 



.(27T) 4 (27r) 4 A(A;)A(g)A(A; + g) 



log.div. 



m 2 ln(m 



(47T) 



: [1 - (4 7 r) 2 a 1 (0)] (76) 



which reproduces the value for d given in Eq. (69) when using a± = —0.030345... [2]. By similar 
considerations we obtain for the logarithmic term of the tensor integral 



a I k d q sin (k^ + q^) 



(2tt)* (2tt)* A(k)A(q)A(k + q y 



= ln(m 2 



log.div. 



(4tt) z 



I + (4^(0) 



(77) 



which again is in agreement with go from Eq. (69). The additional relations derived here are not 
restricted to the leading order in the small mass expansion and one can obtain similar relations 
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for higher order logarithmic pieces. In the remainder of this section, we will briefly illustrate 
the general strategy in relating higher order logarithmic terms of two-loop diagrams to one-loop 
tadpole integrals. For brevity, we will restrict ourselves to ^2 x G(x) 3 , but the generalization to 
tensor integrals is straightforward. We decompose J2 x G(x) 3 as 



$>(x) 3 = G(0) 3 + ]T 

x x^O 



G{xf - F(x) 



(78) 



where the function F(x) is chosen in such a way that the difference in the sum can be trivially 
expanded in the small mass up to a given order. The small mass expansion of F(x) reads 



3 oo 
j=0 i=0 



(79) 



with analytically known expansion coefficients For the n-th order logarithmic term of the 
difference we get 



E 



G(xf - F(x) 



= m 2n \n{m 2 ) £ £ b^a^a^x) - f^x) 

m 2 "ln(m 2 ) X ^ Q l i+j+k=n _ 1 



On the other hand, one has 



Y,x 2l G(x? = 5 lfi G(Q) 2 + Y, 

x x^O 



x Zl G(xY - H(x) 



+ 5>(s), 

x^O 



(80) 



(81) 



where again the function H(x) is chosen in such a way that the difference in the sum can be 
expanded in the small mass up to a given order. If the small mass expansion of H(x) is given 
by 



2 oo 



H(x) = s £ j s £ j h { p(x)m 2i [m 2 \n(m 2 )] j , 

j=0 i=o 

then the l-th order of the difference reads 



E 

x^O 



x 2i G{xf - H{x) 



m 21 J2 Yl x 2i aj(x)a k (x) - hf\x) 

x^O <-j+k=l 



(82) 



(83) 



Since bi(x) is a polynomial in x we can always describe the lattice sum of the logarithmic term 
in Eq. (80) as a linear combination of the terms 



E 



x 2i G{xf - H{x) 



(84) 



and get in this way relations between the logarithmic pieces of Ylx G(x) 3 and terms of the 
type ^2 x x 2%1 . . .x 2ln G(x) 2 . The latter terms can in turn be written as linear combinations of 
tadpoles with tensor structures 



Y,x 2n ■■■X li "G(xf 



[(-D"E 



d \2n 



d \2 



d A k 



(85) 
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such that the logarithmic pieces in J2 x G(x) 3 are expressed in terms of (tensor-)tadpole inte- 
grals. 

6 Conclusions 

In this work, we have extracted the logarithmic and quadratic divergences of the basic sunset 
diagram for a massive scalar field on the lattice and calculated the coefficients of its would-be 
divergences to very high precision. As a first step, the BPHZ fomalism has been applied to the 
sunset integral leading to three integrals without external momentum in the propagators. The 
quadratic and (bi-) logarithmic divergences are then extracted by applying coordinate space 
techniques. 

The crucial observation is that the integrals can be written as products of Green's functions 
and a summation over lattice sites is performed. By subtracting the leading terms of the 
asymptotic form of the Green's function one can expand the products of Green's functions in 
the mass m of the scalar field up to any given order. The coefficients of the logarithmic and 
quadratic divergences are thus expressed in terms of sums over lattice sites. Precise numerical 
knowledge of the Green's function values close to the origin and of its asymptotic large-re 
behavior make an accurate evaluation of these sums possible. Once the external momenta 
have been eliminated from the propagators (e.g. by applying the BPHZ scheme), the method 
proposed can be utilized to calculate any two-loop diagram to very high accuracy. 
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A Asymptotic expansion of G(x) for large x 

In this appendix, we derive the asymptotic expansion of the lattice Green's function G(x) for 
large x. To this aim, it is convenient to introduce the lattice constant a and to study the 
expansion in a with x = ax fixed which is equivalent to the large x = x/a expansion. 
In [6], Paladini and Sexton have studied the asymptotic behavior of the quantity 3 



3 In [6] the results were derived in d dimensions, but here we restrict ourselves to d = 4 dimensions for 
simplicity. 




(A.l) 



where m = a //, and in this appendix we set 




(A.2) 
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In the following, we present an alternative and much easier derivation of their results. The 
small a expansion of the integrand is given by 

1 1 2 a k 4 4 /a(a + l) (A; 4 ) 2 a k 6 

+ a 2 — — , nx ,,, +a 41 



(^2 + /i 2) a (k 2 + /i 2 ) a 12 (k 2 + /i 2 ) a+1 V. 288 (k 2 + jj 2 ) a + 2 360 (k 2 + /i 2 ) a+1 



f " 1 ' 10368 (k 2 + /i 2 )«+ 3 4320 (A; 2 + fx 2 ) a + 2 + 20160 (A; 2 + /i 2 ) Q + 1 
+ 0(a 8 ), (A.3) 

such that the original integral reads 

F a (x) ~ + a 2 Fi 2 )(a;) + a 4 ^(i) + a 6 ^ 6 ^) + 0{a l2 ~ 2a ) (A.4) 

with 

' w/a d 4 k ,^ i 

r/a 

and 



a(a + l)(a + 2) (A; 4 ) 3 a(a + 1) A; 4 A; 6 a 



F ^^C^whr (A - 5) 



*?><*) = £!^(tf )»*£»,(,) + ^fffitt, (A.6) 

*?><«) = a(a+ 10 368 +2> (9|)3F ^W + 2 wWft(') + 20^+' 

The derivative operators in Eqs. (A.6) are defined by d~ = J2t=i^ n /^u- The crucial obser- 
vation is that the extension of the integral boundaries in Eq. (A. 5) to infinity does not affect 
the large- a; behavior of the integral but merely changes its behavior near the origin. For our 
purposes we can thus set for x — > oo 

^^L^whr'^T^^^^' (A - 7) 

where \x\ = y/x 2 . After evaluating the derivatives in Eq. (A.6) we arrive at the results of [6]. 
The asymptotic Green's function G as (x) is immediately obtained from (A.6) by setting a — 1: 

G as (x) * G£)(x) + G<$(x) + G^{x) + G<8(x) + 0(a w ) (A.8) 
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with 

g2>(*) = 

G$(x) = 



1 



An 2 x 2 
1 



m|a;|i ; ri(m|a;|), 



967T 2 (X 2 ) 2 

1 



23040vr 2 (x 2 ) 3 



12(m\x\) K 2 (m\x\) — 6(m|x|) 3 -ft'3(m|a;|) + (m\x\ 



(m\x |) - 2280(m|a;|) 4 K4 



x 



2\2 



KAm\x\ 



+ 180 + 840- 



x 



+ 5(m|x| 



x 



[x 

4\2 



2\2 



(m\x\) 5 K 5 (m\x\) 



60- 



x 



X 2 ) 2 



72- 



x 



x 2 f 



(m\x\) 6 K G (m\x\) 



[x 



2\4 



K 7 {m\x\ 



1 



116121607r 2 (a; 2 ) 4 



ll79360(m\x\) 4 K 4 (m\x\) - 1935360(m\x\) 5 K 5 (m\x\) 



+ (241920 + 1270080-^r) 

v [x z y J 



'm\x\) 6 K e (m\x 



D-( 



7560 + 138600 



x 



(x 2 y 



+ 260064 (^) (m\x\) 7 K 7 (m\x\) + (3780^ + 16380^ + 9072 



x 



X 2 ) 2 
4\2 



(x 2 y 



(x 2 f 



+ lQ200j^y m \x\) 8 K 8 (m\x\) - ^630^ + 1512^) (m|rr|) 9 ^ 9 (m|a;|) 



(x 2 y 



(x 2 y 



+ 35^(m\x\) w K w (m\x\)) 



(A.9) 



The leading asymptotic piece in the small a expansion, G^J , is (up to a prefactor of 1/a 2 ) 
equivalent to the continuum Green's function, while the remaining terms G^\ i > 1, represent 
corrections of order a 2t with respect to G$ and vanish in the continuum limit a — > 0. 
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B Coefficients k and A 

In this appendix, we present the coefficients k and A which appear in Eqs. (29), (35) and (62). 
The coefficients k in Eq. (29) are given in d dimensions by 



«2 



21-3^-6+3^2- f r ( 4 _ 3d) gin (3d.) 

r(2-f) 3 sin(f) ; 
2l -3d m -6+3d n 2-M ( _ 3r (5 _ M) r ^ 2 _ |) + r (4 - ^) T (3 - §)) sin 

r(2-f) 3 r(3-f) sin(f) ; 

2-3^-6+3^2-f ( _ 6r (5 _ M) r (2 - f ) + T (4 - f ) T (3 - f)) sin 



Ko = ^ H K -^- L , B.l 

|) 3 sin(f) 3 

= " ~~(~ 3r ( 5 ~ ~2> l \* - 2) Z t V" ~ ~*> t ^ ~ Ml 0111 (b 2) 



r(2-f) 3 r(3-f)sin(f) 3 
3 ■ 2-^-^ 2 -f r (6 - f ) (r (3 - I ) 2 + 2F (2 - | ) r (4 - j )) sin (g) 



«3 



r(2-|rr(3-rr(4-|)sm(fr 

2 -3<i m -6+3d 7r 2-M (( _ 4 + rf)r (4 _ 3d) + lgr (5 _ 3d)) gin (3|r) 

3(-4 + rf)r(2-f) 3 sin(f) 3 
2 4 ~ M (-22 + 3rf)m- 6+M 7r 2 -¥r (7 - f ) sin (^) 
-8 + d)(-6 + d)(-4 + rf) 2 r (2 - f) 3 sin (f ) 3 



• . . ... . - ,uzjj- (b.4) 



3 • 2 5 - 4d m- 6+3d 7rl-f cot (f ) T (| - f ) 

K4 - ' (R5) 

3 • 2 2 - M m- 6+3 V-¥ cot (f ) ((-4 + d)dT (3 - d) + 4rfr (4 - d) + 2r (5 - d)) 

^ ~ (-2 + d)dY (3 - f) ' (B ' 6) 
3 • 2 2 - M m- 6+M 7r 1 -¥ cot (f ) ((4 - d)d(2 + d)Y (3 - d) + 8(5 + d)r (5 - d))r (4 - f ) 

Kq — k 

{-Q + d)d{-A + d 2 )T (3-|) 2 
3 • 2 4 ~ M (2 + d)m~ 6+3 V-¥ cot (f ) r (5 - d) 

+ (-6 + rf)(-4 + rf 2 )r(3-f) ' (R7) 

_ 3 ■ 2 1 - M m- 6+3 V-¥ 

" 7 " r(|) 2 sin(f) 2 ' 

3 • 2 3 - M m" 6+M 7r 1 -f r (3 - \ ) 

" 8 " (-2 + rf)r(l + |)sin(f) • (R9) 

The coefficients A in Eq. (62) are given by 

A = 0, (B.10) 

Z lib 7T . ... ....... , , , , . , , 

Ai = * ^ * ^ ^ ^ B.ll 

" d\ 3 T fo d\ „• /d7r\ 3 



A 9 = 



T>(5-f) (dT(2-f)+2r (3-f)) sin 

r(2-f) 3 r(3-f) sin(f) ; 

2 -i-3d m -8+3d 7r 2-f r ( 5 _ gd) (rfr (2 - D + 2 r (3 - | )) sin 

d A 3 "P / o d\ I dir\ 3 



r(2-f) d r(3-f) sin(f) 



20 



2-3^-8+3^2-f r ^ _ (dr (2 - |) + 3r (3 - I)) sin (^) 

r(2-|) 2 r(3-|) 2 sin(f) 3 

2 -l-3d dm -8+3d 7r 2-¥ r ( 6 _ M) sin (3Jr) 

r(2-f) 2 r(4-f)sin(f) 3 ' 

2-2-^-8+3^2-^ + 2 cos (d7r))r (5 - f ) (dT (2 - f ) + 2r (3-0) 
r(2-f) 3 r(3-|)sin(f) 2 
2-3^-8+3^2-f (1 + 2 cos (d7r))r (6 - f ) (dT (2 - \ ) + 3r (3 - § )) 

r(2-f) 2 r(3-|) 2 sin(f) 2 

2-1-3^-8+3^2-^ + 2 cos ( d7r )) r ( 7 _ M) (rfr ( 2 - |) + 6r (3 - \ )) 

r(2-f)r(3-f) 3 sin(f) 2 

2-l-3d dm -8+3ci 7r 2-f (! + 2cOs(d7T))r (6 - f ) 

r(2-f) 2 r(4-f)sin(f) 2 

3 • 2' 1 ~ 3d m- s+3d n 2 -'^(l + 2 cos (d7r))r (7 - f ) (dT (2 - §) + T (3 - §)) 

r(2-f) 2 r(3-f)r(4-f)sin(f) 2 

2- 2 - 3d dm' 8+3d 7r 2 ~ f (1 + 2 cos (dn))r (7 - f ) 
r(2-f) 2 r(5-|)sin(f) 2 ' 

-2 2 -^m- 8+M ^-¥ cot r Q - , 

2-i-3d(_ 4 + d)dm- 8+3 V~¥ cot (f ) r (3 - d) 

r(H) 

g-^-s+s^i-M cQt ^ (2 d 2 r ( 4 _ d ) + ( 6 + d ) T (5 - d)) 

(-2 + d)r(3-f) ' 
2 -i-3d(_ 6 + d ) d {2 + d ) m -a+3d n i-3* cot (<|) r(5 _ rf) r(3_|) 

(-2 + d)r(4-f) 2 

2 -3-3d(_ 6 + d )(_ 4 + ^^-8+3^1-f mt ^ p ( 3 _ rf ) 

r(4-|) 

2- M (-6 + d)(-96 + d(-12 + d(3 + d)))m- 8+3 V-ir cot (f ) r (5 - d) 
(-4 + d)(-2 + d)dr (4-f) ' 

2 -3d dm -8+3d 7r 2- f 

r(-i + f)r(f)sin(f) 2 ' 

2 -3d(_ 2 + (_ 4 + d)d) m -8+3d7T 2 -T 

r(f) 2 sm(f) 2 ' 
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